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Ground  State  Resonance  Structure  Calculated  by  Density  Functional  Theory 
For  Estimating  the  Dielectric  Response  of  Some  Typical  High  Explosives 


Introduction 

A  significant  aspect  of  using  response  spectra  calculated  by  density  functional  theory,  DFT,  for 
the  direct  construction  of  permittivity  functions  is  that  it  adopts  the  perspective  of  computational 
physics,  according  to  which  a  numerical  simulation  represents  another  source  of  “experimental”  data. 
This  perspective  is  significant  in  that  a  general  procedure  may  be  developed  for  construction  of 
permittivity  functions  using  DFT  calculations  as  a  quantitative  initial  estimate  of  spectral  response 
features  for  subsequent  adjustment  with  respect  to  additional  information  such  as  experimental 
measurements  and  other  types  of  theory  based  calculations.  That  is  to  say,  for  the  purpose  of  simulating 
many  electromagnetic  response  characteristics  of  materials,  DFT  is  sufficiently  mature  for  the  purpose 
of  generating  data  complementing,  as  well  as  superseding,  experimental  measurements. 

In  the  case  of  THz  excitation  of  materials,  the  procedure  of  using  response  spectra  calculated 
using  DFT  for  the  direct  construction  of  permittivity  functions  is  well  posed  owing  to  the  physical 
characteristic  of  THz  excitation.  In  particular,  it  is  important  to  note  that  the  procedure  for  constructing 
a  permittivity  function  using  response  spectra  calculated  using  DFT  is  physically  consistent  with  the 
characteristically  linear  response  associated  with  THz  excitation  of  molecules.  Accordingly,  one 
observes  a  correlation  between  the  advantages  of  using  THz  excitation  for  detection  of  IEDs  (and 
ambient  materials)  and  those  for  its  numerical  simulation  based  on  DFT.  Specifically,  THz  excitation  is 
associated  with  frequencies  that  are  characteristically  perturbative  to  molecular  states,  in  contrast  to 
frequencies  that  can  induce  appreciable  electronic  state  transitions.  Of  course,  the  practical  aspect  of 
the  perturbative  character  of  THz  excitation  for  detection  is  that  detection  methodologies  can  be 
developed  which  do  not  damage  materials  under  examination.  The  perturbative  character  of  THz 
excitation  with  respect  to  molecular  states  has  significant  implications  with  respect  to  its  numerical 
simulation  based  on  DFT.  It  follows  then  that,  owing  to  the  perturbative  character  of  THz  excitation, 
which  is  characteristically  linear,  one  is  able  to  make  a  direct  association  between  local  oscillations 
about  ground-state  minima  of  a  given  molecule  and  THz  excitation  spectra. 

In  what  follows,  calculations  are  presented  of  ground  state  resonance  structure  associated  with 
the  high  explosives  /?-HMX,  PETN,  RDX,  TNT1  and  TNT2  using  DFT.  This  resonant  structure  is  for 
the  construction  of  parameterized  dielectric  response  functions  for  excitation  by  electromagnetic  waves 
at  compatible  frequencies.  For  this  purpose  the  DFT  software  NRFMOF  was  adopted  [4-10]. 

The  organization  of  the  subject  areas  presented  here  are  as  follows.  First,  a  general  review  of 
the  elements  of  DFT  relevant  for  the  calculation  of  absorption  spectra  is  presented.  This  review  focuses 
on  the  specific  numerical  implementation  of  DFT,  which  is  embodied  by  the  NRFMOF  software. 
Second,  a  general  review  is  presented  concerning  the  formal  structure  of  permittivity  functions  in  terms 
of  analytic  function  representations.  An  understanding  of  the  formal  structure  of  permittivity  functions 
in  terms  of  both  physical  consistency  and  causality  is  important  for  post-processing  of  DFT 
calculations  for  the  purpose  of  constructing  permittivity  functions.  Third,  information  concerning  the 
ground  state  resonance  structure  of  the  explosives  /i-HMX,  PETN,  RDX,  TNT1  and  TNT2,  which  is 
obtained  using  DFT,  is  presented  as  a  set  of  case  studies.  This  information  consists  of  the  ground  state 
molecular  geometry  and  response  spectrum  for  an  isolated  molecule.  In  addition,  for  each  of  the 
explosives,  a  prototype  calculation  is  presented  to  demonstrate  the  construction  of  parameterized 
permittivity  functions  using  response  spectra  calculated  using  DFT.  Fourth,  a  discussion  is  presented 
that  elucidates  the  utility  of  the  information  concerning  the  ground  state  resonance  structure  of  the 
explosives  considered.  This  discussion  also  indicates  the  relevance  of  this  information  for  the 
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construction  of  permittivity  functions  for  frequencies  that  exceed  the  THz  regime.  Finally,  a  conclusion 
is  given,  indicating  possible  future  pathways  for  extension  of  the  methodology  presented  and  the 
calculation  of  spectra  for  other  molecular  systems. 


Construction  of  Permittivity  Functions  using  DFT 


Density  Functional  Theory 

The  application  of  density  functional  theory  (DFT)  and  related  methodologies  for  the 
determination  of  electromagnetic  response  characteristics  is  important  for  the  analyses  of  parameter 
sensitivity.  That  is  to  say,  many  characteristics  of  the  electromagnetic  response  of  a  given  material  may 
not  be  detectable,  or  in  general,  not  relevant  for  detection.  Accordingly,  sensitivity  analyses  concerning 
the  electromagnetic  response  of  layered  composite  systems  can  incorporate  the  results  of  simulations 
using  DFT,  and  related  methodologies,  to  provide  realistic  limits  on  detectability  that  are  independent 
of  a  specific  system  design  for  IED  detection.  In  addition,  analysis  of  parameter  sensitivity  based  on 
atomistic  response  characteristics  of  a  given  material,  obtained  by  DFT,  provide  for  an  “optimal”  best 
fit  of  experimental  measurements  for  the  construction  of  permittivity  functions.  It  follows  that  within 
the  context  of  parameter  sensitivity  analysis,  data  obtained  by  means  of  DFT  represents  a  true 
complement  to  data  that  has  been  obtained  by  means  of  experimental  measurements. 

The  NRLMOL  software  can  be  used  to  compute  an  approximation  of  the  IR  absoprtion 
spectrum  of  a  molecule  [4-10].  NRLMOL  uses  density  functional  theory  to  compute  the  ground  state 
electronic  structure  in  the  Born-Oppenheimer  approximation  using  Kohn-Sham  density  functional 
theory  [2].  NRLMOL  uses  a  Gaussian  orbital  basis  to  describe  the  electronic  wavefunctions  and 
density,  with  numerical  integration  that  is  nearly  exact  to  machine  precision.  Lor  a  given  set  of  nuclear 
positions,  the  calculation  directly  gives  the  electronic  charge  density  of  the  molecule,  the  total  energy 
E ,  and  the  forces  on  each  atom 


dE 

dr' 


(1) 


where  r’a  is  the  a  cartesian  component  of  the  position  of  atom  i,  and  F'a  is  the  corresponding  force. 
The  dipole  moment  of  the  molecule  is  easily  computed  from  the  combined  (nuclear  and  electronic) 
charge  density. 

To  compute  the  miminum  energy  atomic  configuration,  NRLMOL  uses  the  conjugate-gradient 
algorithm  [11].  The  vibrational  spectrum  depends  on  the  atomic  mass  matrix  M iqj/3  =SijSa/3mi  where 

mi  is  the  mass  of  atom  i ,  and  the  energy  second  derivative  matrix 
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oraorp 

through  the  eigenproblem 
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jp 

The  generalized  eigenvalues  £n  are  the  squares  of  the  vibrational  frequencies  con  =  2  m:  vn ,  and  the 
eigenvectors  X("J  give  the  corresponding  atomic  displacements: 
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where  Qn  is  the  normal  mode  coordinate.  The  eigenvectors  are  normalized  according  to  the  condition 

YJx':)MmPx'$=s> .  (5) 

ia.jp 


NRLMOL  computes  the  energy  second  derivatives  Eq.(2)  by  finite  differences,  computing  the  forces 
for  displacement  perturbations  of  each  atom  along  each  Cartesian  direction.  The  first  derivatives  of  the 
dipole  moment  with  respect  to  atomic  positions  djHdr'a  are  also  computed  at  the  same  time.  Each 
vibrational  eigenmode  leads  to  one  peak  in  the  absorption  spectrum,  at  a  frequency  equal  to  the  mode's 
eigenfrequency.  It  is  significant  to  note,  however,  that  the  finite-difference  energy  second  derivatives 
represent  an  approximation  of  the  exact  second  derivatives  and  a  correction  that  reduces  the  associated 
error  of  this  approximation  is  obtained  by  directly  recomputing  the  second  derivatives  of  the  energy 
with  respect  to  the  eigenvectors  displacements. 

The  absorption  intensity  corresponding  to  a  particular  eigenmode  for  a  single  molecule  is  given 
by 


*4  | 
1 

2 

djl 

3c 
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where  c  is  the  speed  of  light  in  a  vacuum,  and 
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The  intensity  Eq.(6)  must  then  be  multiplied  by  the  number  density  of  molecules  to  give  an  absorption 
strength.  It  follows  that  the  absorption  spectrum  calculated  by  NRLMOL  is  a  sum  of  delta  functions 
whose  positions  and  magnitudes  correspond  to  the  vibrational  frequencies  and  magnitudes, 
respectively.  In  principle,  however,  these  spectral  components  must  be  broadened  and  shifted  to 
account  for  anharmonic  effects  such  as  finite  mode  lifetimes  and  inter- mode  couplings. 

Remark.  The  ground  state  resonance  modes  calculated  by  NRLMOL,  which  are  commensurate  with 
electromagnetic  wave  excitation  at  THz  frequencies,  follow  a  “frozen  phonon”  type  method  [12].  This 
method  entails  numerically  a  predictor-corrector  procedure. 

Dielectric  Permittivity  Functions 

The  general  approach  of  constructing  permittivity  functions  according  to  the  best  fit  of  available 
data  for  given  material  corresponding  to  many  different  types  of  experimental  measurements  is  not 
unprecedented  and  has  been  typically  the  dominant  approach,  e.g.,  the  permittivity  function  of  water. 
The  general  simulation  framework  presented  here  considers  an  extension  of  this  approach  in  that 
calculations  of  electromagnetic  response  based  on  DLT  are  incorporated  as  data  for  construction  of 
permittivity  functions.  The  inclusion  of  this  type  of  information  is  essential  for  accessing  what  spectral 
response  features  at  the  molecular  level  are  actually  detectable  with  respect  to  a  given  set  of  detection 
parameters.  Accordingly,  permittivity  functions  having  been  constructed  using  DLT  calculations 
provide  a  quantitative  correlation  between  macroscopic  material  response  and  molecular  structure. 
Within  this  context  it  is  not  important  that  the  permittivity  function  be  quantitatively  accurate  for  the 
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purpose  of  being  adopted  as  input  for  system  simulation.  Rather,  it  is  important  that  the  permittivity 
function  be  qualitatively  accurate  in  terms  of  specific  dielectric  response  features  for  the  purpose  of 
sensitivity  analysis,  which  is  relevant  for  the  assessment  of  absolute  detectability  of  different  types  of 
molecular  structure  with  respect  to  a  given  set  of  detection  parameters.  That  is  to  say,  permittivity 
functions  that  have  been  determined  using  DFT  can  provide  a  mechanistic  interpretation  of  material 
response  to  electromagnetic  excitation  that  could  establish  the  reliability  of  a  given  detection 
methodology  for  detection  of  specific  molecular  characteristics.  Within  the  context  of  practical 
application,  permittivity  functions  having  been  constructed  according  to  the  best  fit  of  available  data 
would  be  “correlated”  with  those  obtained  using  DFT  for  proper  interpretation  of  permittivity-function 
features.  Subsequent  to  establishment  of  good  correlation  between  DFT  and  experiment,  DFT 
calculations  can  be  adopted  as  constraints  for  the  purpose  of  constructing  permittivity  functions,  whose 
features  are  consistent  with  molecular  level  response,  for  adjustment  relative  to  specific  sets  of  either 
experimental  data  or  additional  molecular  level  information. 

The  construction  of  permittivity  functions  using  DFT  calculations  involves,  however,  an  aspect 
that  requires  serious  consideration.  This  aspect  concerns  the  fact  that  a  specific  parametric  function 
representation  must  be  adopted.  Accordingly,  any  parametric  representation,  i.e.,  parameterization, 
adopted  for  permittivity-function  construction  must  be  physically  consistent  with  specific  molecular 
response  characteristics,  while  limiting  the  inclusion  of  feature  characteristics  that  tend  to  mask 
response  signatures  that  may  be  potentially  detectable. 

In  principle,  parameterizations  are  of  two  classes.  One  class  consists  of  parameterizations  that 
are  directly  related  to  molecular  response  characteristics.  This  class  of  parameterizations  would  include 
spectral  scaling  and  width  coefficients.  The  other  class  consists  of  parameterizations  that  are  purely 
phenomenological  and  are  structured  for  optimal  and  convenient  best  fits  to  experimental 
measurements.  A  sufficiently  general  parameterization  of  permittivity  functions  is  given  by  Drude- 
Lorentz  approximation  [13] 


N  yz 

£(V)  =  £(V)  +  i£\v)  =£„+Yj  — - ‘t—. - 


(8) 


where  Vll()  is  the  frequency,  v  and  yn  are  the  spectral  scaling  and  width  of  a  resonance  contributing 
to  the  permittivity  function.  The  permittivity  £rx  is  a  constant  since  the  dielectric  response  at  high 
frequencies  is  substantially  detuned  from  the  probe  frequency.  The  real  £r(V)  and  imaginary  £,(V) 
parts  of  the  permittivity  function  can  be  written  separately  as 
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With  respect  to  practical  application,  the  absorption  coefficient  a  and  index  of  refraction  nr ,  given  by 
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respectively,  provide  direct  relationships  between  calculated  quantities  obtained  by  DFT  and 
“conveniently  measurable”  quantities  a  and  nr . 
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Case  Study  1:  Ground  State  Resonance  Structure  of  /J-HMX 


In  this  section  are  presented  two  sets  of  data,  which  are  the  results  of  computational 
experiments  using  DFT,  concerning  the  molecule  /i-HMX.  These  are  the  relaxed  or  equilibrium 
configuration  of  a  single  isolated  molecule  of  /?-HMX  (see  Table  1)  and  ground-state  oscillation 
frequencies  and  IR  intensities  for  this  configuration  that  are  calculated  by  DFT  according  to  the  frozen 
phonon  approximation  (see  Table  2).  A  schematic  representation  of  the  molecular  geometry  of 
/?-HMX  is  shown  in  Fig.(l). 


Figure  1.  Molecular  Geometry  of  /?-HMX. 

Shown  in  Fig. (2)  is  the  IR  intensity  as  a  function  of  frequency  calculated  using  DFT  for  /9-HMX 
according  to  a  frozen  phonon  approximation.  For  the  spectrum  shown  in  Fig. (2),  the  structure  of  each 
resonance  response  is  approximated  essentially  by  that  of  a  delta  function. 


Table  1.  Atomic  positions  of  /?-HMX  (A) 


Atomic 

number 

X 

Y 

Z 

Atomic 

number 

X 

Y 

Z 

6 

-0.24112 

2.38988 

-0.2519 

8 

-0.51577 

2.12311 

-2.87632 

6 

-2.06734 

0.627 

-0.03867 

8 

-2.13048 

0.63296 

-2.67585 

6 

1.43622 

0.44276 

-0.06484 

8 

1.60928 

3.89177 

0.665 

6 

-0.39044 

-1.31848 

0.14565 

8 

3.1995 

2.36492 

0.62737 

7 

-1.12746 

1.4121 

-0.83025 

8 

-2.23149 

-2.81763 

-0.78797 

7 

1.13398 

1.86281 

-0.12373 

8 

-3.82907 

-1.29822 

-0.73707 

7 

0.49328 

-0.33917 

0.72518 

1 

-0.63362 

2.67459 

0.73503 

7 

-1.76538 

-0.79316 

0.01634 

1 

-0.2064 

3.28046 

-0.88724 

7 

0.62584 

-0.31255 

2.13416 

1 

-3.08482 

0.69866 

-0.44191 

7 

-1.27356 

1.39698 

-2.23846 

1 

-2.04603 

1.05589 

0.97174 

7 

2.05135 

2.75994 

0.44628 

1 

0.00366 

-1.60125 

-0.84102 

7 

-2.67908 

-1.68961 

-0.56012 

1 

-0.42275 

-2.21059 

0.77899 

8 

-0.14183 

-1.03039 

2.77 

1 

2.45279 

0.37224 

0.34102 

8 

1.48106 

0.45171 

2.57423 

1 

1.41841 

0.01191 

-1.07452 
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Table  2.  Oscillation  frequencies  and  IR  intensities 


Frequency 

cm'1 

Intensity 
D2/(amu  A2) 

Frequency 

cm'1 

Intensity 
D2/(amu  A2) 

Frequency 

cm'1 

Intensity 
D2/(amu  A2) 

Frequency 

cm'1 

Intensity 
D2/(amu  A2) 

26.4866 

0.0818 

409.5643 

0.0843 

908.96 

7.8913 

1371.9015 

0.3942 

52.0196 

0.1581 

409.0395 

0.0487 

912.5801 

0.0025 

1388.269 

0.0001 

67.639 

0.0136 

574.8624 

0.0003 

1021.8896 

0.0046 

1411.1066 

0.06 

68.5724 

0.0293 

575.824 

0.8315 

1025.4058 

3.2347 

1411.1876 

1.164 

65.4672 

0.0002 

604.3452 

0.4589 

1098.9735 

4.8331 

1424.5004 

0.0019 

85.7122 

123.1131 

607.001 

0.0003 

1140.1849 

0.0001 

1429.6864 

2.2004 

88.7809 

129.381 

625.655 

0.012 

1179.3564 

0.0008 

1559.9806 

0.0025 

124.5233 

0 

628.2413 

0.0001 

1193.6798 

4.572 

1565.337 

12.9842 

135.4807 

0.007 

708.8374 

0 

1215.4229 

1.9073 

1577.0147 

0.07 

155.4847 

0.0004 

728.9038 

0.3096 

1219.0328 

0.0093 

1580.2385 

13.6073 

158.6671 

0.2157 

731.7556 

0.0004 

1233.6494 

0.0007 

2984.8881 

0.0462 

198.6825 

0.5664 

737.5166 

0.4764 

1241.195 

20.3544 

2985.496 

0.2244 

211.5764 

0.0002 

737.5396 

0.0063 

1257.4516 

7.5375 

2999.5538 

0.7213 

262.3964 

0 

755.5396 

0.2272 

1278.7335 

0.0001 

2998.4535 

0.0822 

283.9432 

0 

818.2214 

0.0067 

1288.5524 

0.001 

3053.2286 

0.0779 

327.8502 

0.1775 

819.2641 

0.1317 

1296.8519 

0.1549 

3053.9098 

0.2106 

330.3013 

0.0001 

847.9616 

1.46 

1313.5215 

0.0001 

3063.1195 

0.0182 

359.6587 

0.1071 

856.2495 

0 

1320.4272 

0.2409 

3063.1191 

0.0105 

388.2644 

0.0013 

894.1733 

5.4895 

1340.3144 

0 

388.1113 

0.1284 

895.9141 

0.5279 

1359.9741 

1.4318 

Figure  2.  IR  intensity  as  a  function  of  frequency  calculated  using  DFT  for  /i-HMX  according  to 
frozen  phonon  approximation. 
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Figure  3.  Real  (blue)  and  imaginary  (red)  parts  of  permittivity  function  of  /?-HMX  molecules  with  yn 
=  3  cm'1  and  p  =  2.4  ■  1019  cm-3  for  frequencies  within  THz  range. 
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Case  Study  2:  Ground  State  Resonance  Structure  of  PETN 


In  this  section  are  presented  two  sets  of  data,  which  are  the  results  of  computational 
experiments  using  DFT,  concerning  the  molecule  PETN.  These  are  the  relaxed  or  equilibrium 
configuration  of  a  single  isolated  molecule  of  PETN  (see  Table  3)  and  ground-state  oscillation 
frequencies  and  IR  intensities  for  this  configuration  that  are  calculated  by  DFT  according  to  the  frozen 
phonon  approximation  (see  Table  4).  A  schematic  representation  of  the  molecular  geometry  of  PETN  is 
shown  in  Fig. (4). 


Figure  4.  Molecular  Geometry  of  PETN. 

Shown  in  Fig. (5)  is  the  IR  intensity  as  a  function  of  frequency  calculated  using  DFT  for  PETN 
according  to  a  frozen  phonon  approximation.  For  the  spectrum  shown  in  Fig. (5),  the  structure  of  each 
resonance  response  is  approximated  essentially  by  that  of  a  delta  function. 


Table  3.  Atomic  positions  of  PETN  (A). 


Atomic 

number 

X 

Y 

Z 

Atomic 

number 

X 

Y 

Z 

6 

-0.90787 

0.48907 

0.17354 

8 

1.94808 

-3.14586 

-0.03192 

6 

0.54956 

0.01783 

0.36869 

8 

2.80821 

-1.21442 

0.58616 

6 

-0.96651 

2.01269 

0.41718 

1 

-1.99624 

2.37482 

0.29107 

6 

-1.80732 

-0.27726 

1.16724 

1 

-0.31204 

2.54022 

-0.28938 

6 

-1.40666 

0.20317 

-1.25944 

1 

0.87858 

0.19321 

1.40132 

8 

-0.52463 

2.24139 

1.76964 

1 

1.21509 

0.56308 

-0.3142 

8 

0.58244 

-1.39141 

0.07067 

1 

-1.46281 

-0.10281 

2.19548 

8 

-3.14643 

0.22812 

0.99941 

1 

-1.78105 

-1.35577 

0.96103 

8 

-0.5109 

0.88215 

-2.1613 

1 

-2.43196 

0.57423 

-1.39229 

8 

-5.23155 

0.00318 

1.74575 

1 

-1.38953 

-0.87791 

-1.45237 

8 

-3.67634 

-1.26381 

2.65601 

7 

-4.11793 

-0.42256 

1.90605 

8 

-0.08736 

3.85131 

3.24351 

7 

-0.45302 

3.67996 

2.11044 

8 

-0.75847 

4.4583 

1.23535 

7 

1.93416 

-1.97171 

0.22886 

00  OO 

-0.09204 

-1.81642 

1.22306 

-0.04519 

-4.32045 

-3.80006 

7 

-0.85337 

0.65487 

-3.58248 
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Table  4.  Oscillation  frequencies  and  IR  intensities 


Frequency 

cm'1 

Intensity 
D2/(amu  A2) 

Frequency 

cm'1 

Intensity 
D2/(amu  A2) 

Frequency 

cm'1 

Intensity 
D2/(amu  A2) 

Frequency 

cm'1 

Intensity 
D2/(amu  A2) 

27.6924 

0.0072 

307.6541 

0 

893.5171 

0.3109 

1341.9832 

0.9802 

45.7692 

0.0327 

419.9976 

0.2946 

908.4957 

0.1288 

1345.3347 

1.6302 

75.5706 

0.0023 

420.3111 

0.2982 

907.8127 

0.135 

1358.2592 

0.0004 

65.4508 

0.015 

483.3724 

1.6977 

973.4102 

0.0007 

1450.6195 

0.0268 

82.2402 

0.0111 

555.0496 

0.0005 

995.1452 

1.022 

1451.8016 

0.0914 

85.8255 

0.0045 

577.1765 

0.8504 

995.8256 

1.0268 

1456.0802 

0.397 

52.9438 

0.0048 

579.001 

0.8601 

1024.0042 

0.0008 

1456.6128 

0.3994 

78.9206 

0.0323 

583.1426 

0.2211 

1033.1369 

1.8175 

1690.1603 

1.0193 

77.8851 

0.0334 

609.0535 

0.0006 

1133.8946 

0.0215 

1692.9065 

7.8371 

124.1609 

0.0096 

652.6221 

2.8048 

1157.0829 

0.0091 

1693.0525 

12.4003 

124.381 

0.0072 

651.6133 

2.7891 

1157.2312 

0.0092 

1697.4129 

11.0743 

145.0029 

0.0471 

704.123 

4.6868 

1216.6553 

0.0005 

2982.9778 

0.1752 

146.806 

0.0001 

727.4008 

0.1249 

1234.5454 

0.494 

2985.3071 

0.2 

171.1139 

0.0014 

726.4067 

0.1068 

1235.3559 

0.514 

2984.8943 

0.1837 

192.1037 

0.0366 

727.571 

0.1801 

1254.1184 

6.7697 

2986.1399 

0.0196 

192.0644 

0.0295 

728.1667 

0.8163 

1267.8336 

6.9384 

3035.5243 

0.0363 

196.0122 

0.0006 

787.2791 

7.1124 

1267.4456 

6.9192 

3038.1006 

0.0998 

237.8431 

0.0229 

787.7443 

7.1507 

1280.9494 

0.0544 

3038.8521 

0.1325 

240.7368 

0.014 

788.9399 

12.423 

3039.7682 

0.2253 

242.4609 

0.0059 

802.3902 

0.0015 

1284.8749 

4.0339 

295.694 

0.0012 

825.5549 

0.0012 

1342.235 

1 

Wavenumber  (2tz/X)  (rad/cm) 

Figure  5.  IR  intensity  as  a  function  of  frequency  calculated  using  DFT  for  PETN  according  to  frozen 
phonon  approximation. 
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Figure  6.  Real  (blue)  and  imaginary  (red)  parts  of  permittivity  function  of  PETN  molecules  with  yn 
3  cm  1  and  p  =  2.4  ■  1019  cm  for  frequencies  within  THz  range. 
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Case  Study  3:  Ground  State  Resonance  Structure  of  RDX 


In  this  section  are  presented  two  sets  of  data,  which  are  the  results  of  computational 
experiments  using  DFT,  concerning  the  molecule  RDX.  These  are  the  relaxed  or  equilibrium 
configuration  of  a  single  isolated  molecule  of  RDX  (see  Table  5)  and  ground-state  oscillation 
frequencies  and  IR  intensities  for  this  configuration  that  are  calculated  by  DFT  according  to  the  frozen 
phonon  approximation  (see  Table  6).  A  schematic  representation  of  the  molecular  geometry  of  RDX  is 
shown  in  Fig. (7). 


Figure  7.  Molecular  Geometry  of  RDX. 


Shown  in  Fig. (8)  is  the  IR  intensity  as  a  function  of  frequency  calculated  using  DFT  for  RDX 
according  to  a  frozen  phonon  approximation.  For  the  spectrum  shown  in  Fig. (8),  the  structure  of  each 
resonance  response  is  approximated  essentially  by  that  of  a  delta  function. 


Table  5.  Atomic  positions  of  RDX  (A): 


Atomic 

number 

X 

Y 

Z 

Atomic 

number 

X 

Y 

Z 

7 

1.16994 

-0.81313 

0.09076 

8 

-1.12995 

1.87973 

2.52487 

7 

0.09076 

1.16994 

-0.81313 

8 

3.34946 

-0.19764 

0.06923 

7 

-0.81313 

0.09076 

1.16994 

8 

0.06923 

3.34946 

-0.19764 

7 

2.44674 

-0.67765 

0.7438 

8 

-0.19764 

0.06923 

3.34946 

7 

0.7438 

2.44674 

-0.67765 

1 

-1.50333 

1.88316 

0.33719 

7 

-0.67765 

0.7438 

2.44674 

1 

0.33719 

-1.50333 

1.88316 

6 

-1.09011 

0.92662 

0.00735 

1 

1.88316 

0.33719 

-1.50333 

6 

0.00735 

-1.09011 

0.92662 

1 

-1.82611 

0.39667 

-0.61584 

6 

0.92662 

0.00735 

-1.09011 

1 

-0.61584 

-1.82611 

0.39667 

00  OO 

2.52487 

1.87973 

-1.12995 

2.52487 

1.87973 

-1.12995 

1 

0.39667 

-0.61584 

-1.82611 
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Table  6.  Oscillation  frequencies  and  IR  intensities 


Frequency 

cm'1 

Intensity 
D2/(amu  A2) 

Frequency 

cm"1 

Intensity 
D2/(amu  A2) 

Frequency 

cm"1 

Intensity 
D2/(amu  A2) 

Frequency 

cm"1 

Intensity 
D2/(amu  A2) 

50.4476 

0.0012 

561.592 

0.068 

958.8454 

0.841 

1394.1574 

0.5447 

51.1589 

0.0012 

561.6381 

0.0682 

958.8824 

0.8369 

1395.0201 

0.5394 

77.7345 

0.0128 

573.1913 

0.0001 

1106.4304 

0 

1414.4218 

1.9899 

91.1576 

0.0081 

635.6372 

0.2794 

1185.1197 

0.598 

1572.3865 

0.0002 

125.232 

0.0004 

635.8515 

0.2799 

1208.9255 

2.0878 

1600.7625 

8.9248 

125.1679 

0.0004 

723.8626 

0.0454 

1209.0906 

2.0855 

1600.6959 

8.9225 

215.3785 

0.1185 

728.3136 

0.0363 

1218.0051 

0.0005 

2979.9541 

0.0238 

215.3755 

0.1184 

728.3277 

0.036 

1231.9493 

2.5768 

2980.1722 

0.0204 

290.8681 

0.0001 

752.231 

1.6486 

1231.8716 

2.5783 

2983.6473 

1.0283 

346.3176 

0.0151 

830.6969 

0.908 

1284.8248 

6.3921 

3098.3191 

0.3366 

346.2359 

0.0151 

830.9032 

0.9079 

1304.8515 

0.0003 

3099.4052 

0.3341 

392.0173 

0.0047 

859.1975 

6.4437 

1326.7768 

0.368 

3096.9028 

0.4219 

391.356 

0.0047 

858.8267 

6.4375 

1327.4575 

0.3693 

413.5828 

0.3457 

863.9551 

0.0141 

1344.0359 

0.6348 

431.3056 

0.5958 

886.788 

3.5023 

1344.7081 

0.6313 

Figure  8.  IR  intensity  as  a  function  of  frequency  calculated  using  DFT  for  RDX  according  to  frozen 
phonon  approximation. 
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Figure  9.  Real  (blue)  and  imaginary  (red)  parts  of  permittivity  function  of  RDX  molecules  with  yn  =  3 
cm  1  and  p  =  2.4- 1019  cm  3  for  frequencies  within  THz  range. 
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Case  Study  4:  Ground  State  Resonance  Structure  of  TNT1 


In  this  section  are  presented  two  sets  of  data,  which  are  the  results  of  computational 
experiments  using  DFT,  concerning  the  molecule  TNT1.  These  are  the  relaxed  or  equilibrium 
configuration  of  a  single  isolated  molecule  of  TNT1  (see  Table  7)  and  ground-state  oscillation 
frequencies  and  IR  intensities  for  this  configuration  that  are  calculated  by  DFT  according  to  the  frozen 
phonon  approximation  (see  Table  8).  A  schematic  representation  of  the  molecular  geometry  of  TNT1  is 
shown  in  Fig. (10). 


Figure  10.  Molecular  Geometry  of  TNT  1. 

Shown  in  Fig. (11)  is  the  IR  intensity  as  a  function  of  frequency  calculated  using  DFT  for  TNT1 
according  to  a  frozen  phonon  approximation.  For  the  spectrum  shown  in  Fig.(ll),  the  structure  of  each 
resonance  response  is  approximated  essentially  by  that  of  a  delta  function. 


Table  7.  Atomic  positions  of  TNT1  (A): 


Atomic 

number 

X 

Y 

Z 

Atomic 

number 

X 

Y 

Z 

1 

1.277 

-0.71765 

1.9543 

6 

1.41827 

-0.47815 

-2.75813 

1 

-2.49724 

-0.7602 

-0.12405 

7 

-1.38812 

-0.82281 

2.33281 

1 

0.77359 

-0.73572 

-3.60038 

7 

2.83728 

-0.47406 

-0.14338 

1 

2.29092 

-1.14373 

-2.77796 

7 

-1.56142 

-0.74062 

-2.57029 

1 

1.80338 

0.54246 

-2.88489 

8 

3.37643 

0.34145 

-0.89378 

6 

0.71266 

-0.68704 

1.02471 

8 

3.41117 

-1.19787 

0.67094 

6 

-0.67498 

-0.7428 

1.03277 

8 

-0.69699 

-0.7987 

3.35263 

6 

-1.40976 

-0.72761 

-0.14562 

8 

-2.61686 

-0.90672 

2.29484 

6 

-0.71255 

-0.68459 

-1.3493 

8 

-2.49644 

0.05804 

-2.62488 

6 

6 

0.69021 

1.35505 

-0.58801 

-0.58239 

-1.44921 

-0.20393 

8 

-1.28456 

-1.59556 

-3.41343 
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Table  8.  Oscillation  frequencies  and  IR  intensities 


Frequency 

cm'1 

Intensity 
D2/(amu  A2) 

Frequency 

cm'1 

Intensity 
D2/(amu  A2) 

Frequency 

cm'1 

Intensity 
D2/(amu  A2) 

Frequency 

cm'1 

Intensity 
D2/(amu  A2) 

172.5174 

0.0237 

414.8762 

0.0385 

908.0919 

0.4081 

1449.4697 

0.2159 

197.4056 

0.011 

495.2412 

0.0797 

913.9983 

0.5736 

1487.7755 

0.4935 

135.2506 

0.001 

522.6106 

0.0393 

1057.9824 

0.0914 

1528.3766 

4.1144 

284.1745 

0.0427 

560.2942 

0.0077 

1068.9096 

0.1036 

1524.3186 

6.5037 

361.7698 

0.012 

578.169 

0.0073 

1057.3348 

1.2247 

1532.7054 

3.8635 

427.1822 

0.133 

668.3243 

0.432 

1137.2714 

0.3466 

1585.1704 

1.0497 

172.7674 

0.0458 

702.98 

0.0633 

1177.5714 

0.2562 

1592.6665 

0.7857 

171.747 

0.0908 

709.4691 

1.031 

1190.6768 

0.0231 

2804.9617 

0.0052 

177.736 

0.0114 

733.5913 

0.398 

1301.8607 

8.6558 

2860.6966 

0.0794 

285.6632 

0.0627 

744.7193 

0.1154 

1308.2775 

7.0606 

2930.3248 

0.053 

300.4711 

0.0241 

768.8811 

0.017 

1325.0594 

0.5661 

3154.161 

0.4024 

304.604 

0.0077 

783.0319 

0.2198 

1355.3609 

0.1079 

3160.6014 

0.5166 

337.6068 

0.0188 

808.5178 

0.0437 

1385.5683 

0.2682 

342.8559 

0.0512 

889.7004 

0.8604 

1386.0299 

0.0801 

376.9598 

0.053 

907.8869 

0.4586 

1462.453 

0.0542 

Wavenumber  (1/M  (cm'1) 

Figure  11.  IR  intensity  as  a  function  of  frequency  calculated  using  DFT  for  TNT1  according  to  frozen 
phonon  approximation. 
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Figure  12.  Real  (blue)  and  imaginary  (red)  parts  of  permittivity  function  of  TNT1  molecules  with  yn 
3  cm  1  and  p  -  2.4  ■  1019  cm  for  frequencies  within  THz  range. 
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Case  Study  5:  Ground  State  Resonance  Structure  of  TNT2 


In  this  section  are  presented  two  sets  of  data,  which  are  the  results  of  computational 
experiments  using  DFT,  concerning  the  molecule  TNT2.  These  are  the  relaxed  or  equilibrium 
configuration  of  a  single  isolated  molecule  of  TNT2  (see  Table  9)  and  ground-state  oscillation 
frequencies  and  IR  intensities  for  this  configuration  that  are  calculated  by  DFT  according  to  the  frozen 
phonon  approximation  (see  Table  10).  A  schematic  representation  of  the  molecular  geometry  of  TNT2 
is  shown  in  Fig.(13). 


Figure  13.  Molecular  Geometry  of  TNT2. 

Shown  in  Fig.(14)  is  the  IR  intensity  as  a  function  of  frequency  calculated  using  DFT  for  TNT2 
according  to  a  frozen  phonon  approximation.  For  the  spectrum  shown  in  Fig. (14),  the  structure  of  each 
resonance  response  is  approximated  essentially  by  that  of  a  delta  function. 


Table  9.  Atomic  positions  of  TNT2  (A): 


Atomic 

number 

X 

Y 

Z 

Atomic 

number 

X 

Y 

Z 

1 

1.277 

-0.71765 

1.9543 

6 

1.41827 

-0.47815 

-2.75813 

1 

-2.49724 

-0.7602 

-0.12405 

7 

-1.38812 

-0.82281 

2.33281 

1 

0.77359 

-0.73572 

-3.60038 

7 

2.83728 

-0.47406 

-0.14338 

1 

2.29092 

-1.14373 

-2.77796 

7 

-1.56142 

-0.74062 

-2.57029 

1 

1.80338 

0.54246 

-2.88489 

8 

3.37643 

0.34145 

-0.89378 

6 

0.71266 

-0.68704 

1.02471 

8 

3.41117 

-1.19787 

0.67094 

6 

-0.67498 

-0.7428 

1.03277 

8 

-0.69699 

-0.7987 

3.35263 

6 

-1.40976 

-0.72761 

-0.14562 

8 

-2.61686 

-0.90672 

2.29484 

6 

-0.71255 

-0.68459 

-1.3493 

8 

-2.49644 

0.05804 

-2.62488 

6 

6 

0.69021 

1.35505 

-0.58801 

-0.58239 

-1.44921 

-0.20393 

8 

-1.28456 

-1.59556 

-3.41343 
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Table  10.  Oscillation  frequencies  and  IR  intensities 


Frequency 

cm'1 

Intensity 
D2/(amu  A2) 

Frequency 

cm'1 

Intensity 
D2/(amu  A2) 

Frequency 

cm"1 

Intensity 
D2/(amu  A2) 

Frequency 

cm’1 

Intensity 
D2/(amu  A2) 

172.5174 

0.0237 

414.8762 

0.0385 

908.0919 

0.4081 

1449.4697 

0.2159 

197.4056 

0.011 

495.2412 

0.0797 

913.9983 

0.5736 

1487.7755 

0.4935 

135.2506 

0.001 

522.6106 

0.0393 

1057.9824 

0.0914 

1528.3766 

4.1144 

284.1745 

0.0427 

560.2942 

0.0077 

1068.9096 

0.1036 

1524.3186 

6.5037 

361.7698 

0.012 

578.169 

0.0073 

1057.3348 

1.2247 

1532.7054 

3.8635 

427.1822 

0.133 

668.3243 

0.432 

1137.2714 

0.3466 

1585.1704 

1.0497 

172.7674 

0.0458 

702.98 

0.0633 

1177.5714 

0.2562 

1592.6665 

0.7857 

171.747 

0.0908 

709.4691 

1.031 

1190.6768 

0.0231 

2804.9617 

0.0052 

177.736 

0.0114 

733.5913 

0.398 

1301.8607 

8.6558 

2860.6966 

0.0794 

285.6632 

0.0627 

744.7193 

0.1154 

1308.2775 

7.0606 

2930.3248 

0.053 

300.4711 

0.0241 

768.8811 

0.017 

1325.0594 

0.5661 

3154.161 

0.4024 

304.604 

0.0077 

783.0319 

0.2198 

1355.3609 

0.1079 

3160.6014 

0.5166 

337.6068 

0.0188 

808.5178 

0.0437 

1385.5683 

0.2682 

342.8559 

0.0512 

889.7004 

0.8604 

1386.0299 

0.0801 

376.9598 

0.053 

907.8869 

0.4586 

1462.453 

0.0542 

Figure  14.  IR  intensity  as  a  function  of  frequency  calculated  using  DFT  for  TNT2  according  to  frozen 
phonon  approximation. 
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0.016 


Figure  15.  Real  (blue)  and  imaginary  (red)  parts  of  permittivity  function  of  TNT2  molecules  with  yn 
3  cm  1  and  p  =  2.4  ■  1019  cm  for  frequencies  within  THz  range. 
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Discussion 


The  DFT  calculated  absorption  spectra  given  in  tables  2,  4,  6  and  8  provide  two  types  of 
information  for  general  analysis  of  dielectric  response.  These  are  the  denumeration  of  ground  state 
resonance  modes  and  estimates  of  molecular  level  dielectric  response  structure.  The  construction  of 
permittivity  functions  using  the  DFT  calculated  absorption  spectra  follows  the  same  procedure  as  that 
applied  for  the  construction  of  permittivity  functions  using  experimentally  measured  absorption 
spectra,  but  with  the  addition  of  certain  constraint  conditions.  Accordingly,  construction  of  permittivity 
functions  using  either  DFT  or  experimentally  measured  absorption  spectra  requires  parameterizations 
that  are  in  terms  of  physically  consistent  analytic  function  representations  such  as  the  Drude-Lorentz 
model.  Although  the  formal  structure  of  permittivity  functions  constructed  using  DFT  and  experimental 
measurements  are  the  same,  their  interpretation  with  respect  to  parameterization  is  different  for  each 
case. 

The  construction  of  permittivity  functions  using  experimental  measurements,  an  established 
methodology,  defines  an  inverse  problem  where  resonant  locations,  peaks  and  widths,  as  well  as  the 
number  of  resonances,  are  assumed  adjustable.  Following  this  approach,  it  follows  that  many  of  the 
detailed  characteristics  of  resonance  structure  are  smoothed  or  averaged.  In  addition,  measurement 
artifacts  associated  with  sample  preparation  and  detector  configuration  can  in  principle  introduce 
errors.  One  advantage  of  permittivity  functions  constructed  using  experimental  measurements, 
however,  is  that  many  aspects  of  dielectric  response  on  the  macroscale  that  are  associated  with 
multiscale  averaging  and  molecule-lattice  coupling  are  taken  into  account  inherently.  Accordingly,  the 
disadvantage  of  this  approach  is  that  the  nature  of  any  multiscale  averaging  and  resonant  structure, 
contributing  to  dielectric  response  on  the  macroscopic  level,  may  not  be  understood.  This  lack  of 
quantitative  understanding  can  in  principle  inhibit  the  development  of  pump-probe  type  methodologies 
for  selective  excitation  of  molecular  modes,  which  are  for  the  purpose  of  enhanced  signature  detection 
or  modulation. 

The  construction  of  permittivity  functions  using  DFT  calculations,  the  methodology  whose 
development  is  considered  here,  defines  a  direct  problem  approach  where  dielectric  response  is 
estimated  within  the  bounds  of  relatively  well-defined  adjustable  parameters.  Following  this  approach, 
a  peimittivity  function  is  constructed  using  the  DFT  calculated  absorption  spectra,  e.g.,  tables  2,  4,  6,  8 
and  10,  under  the  condition  that  the  calculated  resonance  locations  are  fixed,  while  resonance  widths 
and  number  densities  are  assumed  adjustable,  e.g.,  Figs.  (3),  (6),  (9),  (12)  and  (15).  Better  interpretation 
of  dielectric  response  of  explosives  on  a  macroscale  can  be  achieved  through  correlation  of  resonance 
structure  that  is  experimentally  observed  and  calculated  by  DFT.  In  principle,  correlation  of  resonance 
structure  would  include  the  quantitative  analysis  of  changes  in  signature  features  associated  with  the 
transition  of  the  system  from  that  of  a  low-density  system  of  uncoupled  molecule  to  that  of  a  bulk 
lattice. 

With  respect  to  more  extensive  DFT  calculations  concerning  the  ground-state  absorption  spectra 
of  a  bulk  lattice  or  spectra  corresponding  to  electronic  state  transitions,  it  is  important  to  note  that  the 
atomic  positions  of  the  relaxed  or  equilibrium  configuration  of  a  single  isolated  molecule,  e.g.,  tables  1, 
3,  5,  7  and  9,  provide  a  convenient  starting  point.  Calculation  of  the  dielectric  response  of  a  bulk  lattice 
would  entail,  in  principle,  the  construction  of  a  super  cell  consisting  of  molecules  whose  initial 
positions  are  those  determined  by  DFT  for  isolated  systems.  Additional  constraints  on  this  super  cell 
could  be  based  on  crystallographic  information  concerning  bulk  density  or  lattice  spacing.  Calculation 
of  the  dielectric  response  associated  with  electronic  state  transitions  would  entail  the  application  of 
methods  based  on  perturbation  theory.  In  principle,  for  these  methods  most  of  the  computational  effort 
is  expended  for  determination  of  the  ground  state,  with  respect  to  which  all  excited  states  are 
determined  self  consistently.  These  methods  typically  would  be  based  on  time-dependent  density 
functional  theory  (TDDFT). 
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Conclusion 


The  calculations  of  ground  state  resonance  structure  associated  with  the  high  explosives 
/?-HMX,  PETN,  RDX,  TNT1  and  TNT2  using  DFT  are  meant  to  serve  as  reasonable  estimates  of 
molecular  level  response  characteristics,  providing  interpretation  of  dielectric  response  features,  for 
subsequent  adjustment  relative  to  experimental  measurements  and  molecular  structure  theory. 
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Appendix  1 


General  Framework  for  Numerical  Simulation  of  IED  Detection  and 
Remote  Activation  Scenarios 

The  set  of  parameterized  permittivity  functions  given  above  represent  contributions  to  a 
component  model  of  a  general  framework  for  numerical  simulation  of  IED  detection  and  remote 
activation  scenarios,  whose  initial  constmction  was  described  in  reference  14  (see  Fig.(Al)).  In 
addition,  the  discussions  concerning  NRLMOL  given  above,  as  well  as  previously  [14],  provide 
information  concerning  the  practical  application  of  NRLMOL  for  construction  of  dielectric  response 
functions.  Accordingly,  included  in  this  appendix  are  notes  describing  modifications  and  errata 
associated  with  the  continuous  evolution  of  this  framework  for  purposes  of  practical  analysis  and 
simulation.  In  particular,  included  in  this  appendix  are  corrections  to  the  S-matrix  computer  program 
presented  in  reference  14. 


I -  i 


Figure  Al.  General  framework  for  numerical  simulation  of  IED  detection. 
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Discussion  Concerning  General  Framework  for  Numerical  Simulation  of  IED  detection. 

The  general  framework  for  numerical  simulation  of  IED  detection  and  remote  activation 
scenarios  presented  in  reference  14  places  emphasis  on  the  construction  of  permittivity  functions  using 
DFT  calculated  spectra,  their  parameterization  and  parameter  adjustment  with  respect  to  additional 
information,  which  are  for  input  to  an  S-matrix  representation  of  a  layered  and  composite  system.  It  is 
significant  to  note,  however,  that  from  a  practical  perspective,  the  construction  of  permittivity  functions 
using  DFT  represents  one  of  three  approaches  to  inputting  information  concerning  dielectric  response 
of  materials  to  a  general  framework  based  on  an  S-matrix  representation.  These  other  approaches  are: 
construction  of  permittivity  functions  by  inverse  analysis  of  experimentally  measured  spectra;  and  the 
use  of  permittivity  functions  that  have  been  previously  reported  and  documented.  Accordingly,  a 
schematic  of  the  general  framework  for  numerical  simulation  of  IED  detection  is  shown  in  Fig.(Al), 
which  describes  the  interrelation  between  three  approaches  for  inputting  information  concerning 
dielectric  response  of  materials  to  the  S-matrix  representation.  Referring  to  Fig.(Al),  it  is  to  be  noted 
that  construction  of  permittivity  functions  using  DFT  calculated  spectra  and  experimentally  measured 
spectra  follow  approaches  that  are  similar,  and  therefore  represent  a  single  component  of  the  general 
framework,  which  is  shown  as  contained  within  dashed  lines.  Referring  again  to  Fig.(Al),  it  is  seen  that 
from  the  perspective  of  practical  analysis,  a  data  base  consisting  of  previously  reported  and 
documented  permittivity  functions  represents  a  component  of  the  general  framework  that  extends  its 
capability  with  respect  to  simulation  of  complex  systems.  In  particular,  the  response  signatures  of 
permittivity  functions  constructed  using  DFT  calculated  spectra  and  experimentally  measured  spectra 
can  be  examined  as  they  occur  in  combination  with  materials  whose  dielectric  response  has  already 
been  well  documented  quantitatively,  e.g.,  metals  [15]  and  water  [16,17]. 
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Itemized  corrections  of  figures  presented  in  reference  14. 

In  Fig. (3)  of  reference  14,  the* -axis  scale  should  be  multiplied  by  2k  for  correct  values  of  the 
wavenumber.  In  addition,  this  figure  does  not  shows  the  resonances  bellow  300  cin(~  50x2#  cm'1). 
For  frequencies  lower  than  300  cm'1,  the  permittivity  function  of  /3  -HMX  is  shown  in  Fig.(A2). 


Figure  A2.  Permittivity  function  of  /J-HMX  for  frequencies  within  THz  range. 
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In  Fig.(4)  of  reference  14,  the  x  -axis  scale  should  be  multiplied  by  2 n  for  correct  values  of  the 
wavenumber.  For  frequencies  lower  than  300  cm_1(  ~  50x2#  cm'1),  the  reflectivity  function  Rs  is 
shown  in  Fig. (A3)  for  an  optically  dense  layer  of  /5-HMX  molecules  with  =  1  (upper)  and  =  1.2 
(lower). 


Figure  A3.  Reflectivity  functions  for  a  layer  of  /?-HMX  on  a  gold  substrate. 
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In  Figs.  (7)  and  (8)  of  reference  14,  the  x  -axis  scale  should  be  multiplied  by  2 n  for  correct 
values  of  the  wavenumber.  For  frequencies  lower  than  300  cm  ’(=  50x  2 n  cm"1),  the  absorption 
coefficient  and  index  of  refraction  are  shown  in  the  Fig.(A4)  for  y  =  3  cm  and  N  =  4  cm  *. 


Figure  A4.  Absorption  coefficient  and  index  of  refraction  for  /i-HMX  calculated  by  DFT  for  THz 

12 

range  of  frequencies,  where  adjustable  parameters  y  =  3  cm  and  N  =  4  cm"'. 
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Table  Al.  Modifications  to  S-matrix  computer  program  presented  in  reference  14. 


Page  # 

Old 

New 

16 

if  (ISCAN  .eq.  0)  then  !  angle  scan 

write(2,*) '  ANGLE  Rp  Rs  ' 

do  I  =  LNMESHA+NBINA-1 

if  (ISCAN  .eq.  0)  then  !  angle  scan 

write(2,*) '  ANGLE  Rp  Rs  ' 

do  I  =  LNMESHA 

16 

write(3,*) 1  ANGLE  Rp  Rs  ’ 

do  I  =  1,NMESHA 

write(3,*) '  ANGLE  Rp  Rs  ' 

do  I  =  1,NMESHA-NBINA2 

16 

do  IB  =  1,NBINA 

reft  b=ref  1  b+ref  1  (I-NBIN  A2+IB  - 1 ) 
ref2b=ref2b+ref2(I-NBINA2+IB-l ) 
end  do 

do  IB  =  1,NBINA 

ref  1  b=ref  lb+ref  1  (I-NBIN  A2+IB) 
ref2b=ref2b+ref2(I-NBINA2+IB) 
end  do 

17 

else  !  frequency  scan 

write(2,*) '  wnum  frequency  Rp  Rs  ' 

do  J  =  LNMESHF+NBINF-1 

else  !  frequency  scan 

write(2,*) '  wnum  frequency  Rp  Rs  ' 

do  J  =  1,NMESHF 

17 

write(3,*) '  wnum  frequency  Rp  Rs  ' 

do  J  =  1,NMESHF 
reflb=0.0d0 
ref2b=0.0d0 
do  JB  =  LNBINF 

ref  1  b=ref  1  b+ref  1  (J-NBINF2+ JB- 1 ) 
ref2b=ref2b+ref2(J-NBINF2+JB-l) 
end  do 

write(3,*) '  wnum  frequency  Rp  Rs  ' 

do  J  =  NBINF2.NMESHF-NBINF2 
reflb=0.0d0 
ref2b=0.0d0 
do  JB  =  LNBINF 
ref  lb=ref  1  b+ref  1  (J-NBINF2+JB) 
ref2b=ref2b+ref2(J-NBINF2+JB) 
end  do 
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